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1.  Introduction 


Electromagnetic-oriented  detection  is  expected  to  be  an  integral  component  of  the  U.S.  Army’s 
multi-sensor-based  strategy  for  force  protection  and  route  clearance.  Over  the  last  decade,  a 
myriad  of  research  avenues  have  been  directed  toward  developing  ground-based,  ultra- wideband 
(UWB)  radars  with  the  capability  to  simultaneously  penetrate  the  ground,  and  image  concealed 
landmines  and  improvised  explosive  devices  (IEDs)  from  a  standoff  distance  (Kositsky  et  al., 
2002;  Bradley  et  al.,  2003;  Wang  et  al.,  2005;  Ressler  et  al.,  2007;  Wang  et  al.,  2008).  Despite  the 
fact  that  experimentation  over  the  last  several  years  has  demonstrated  the  potential  of  the 
technology  against  certain  classes  of  targets,  the  amorphous  nature  of  the  encountered  threats  has 
complicated  the  development  and  fielding  of  the  technology.  One  promising  system — developed 
at  the  U.S.  Army  Research  Laboratory  (ARL) — is  the  low-frequency,  UWB  synchronous  impulse 
reconstruction  (SIRE),  imaging  ground-penetrating  radar  (GPR)  (figure  1),  which  is  a 
vehicle-based  platform  operating  at  the  nominal  frequency  band  of  500-2500  MHz  with  the 
forward-looking  coverage  angle  spanning  approximately  5°-20°  with  respect  to  the  horizon.  As 
opposed  to  conventional  downward-looking  GPR,  a  forward-looking  GPR  such  as  the  SIRE 
system  offers  the  advantage  of  having  a  reduced  ground-bounce  return  due  to  the  use  of  a  sensing 
geometry  at  an  off-normal — or  near-grazing — angle.  One  of  the  disadvantages  of  the 
forward-looking  geometry  is  realized  by  first  noting  that  the  overall  backscattered  signal  power 
response  of  a  target,  either  on-surface  or  buried,  decreases  as  a  function  of  the  illumination- angle 
approximately  proportional  to  the  form  of  \(l+R(6uiu))\4,  where  R{0uh )  is  the  reflection  coefficient 
of  the  ground-surface.  (Here,  it  is  assumed  the  predominant  interactions  are  from  plane-wave 
far-field  effects — with  negligible  surface-wave  contributions).  This  behavior  of  the  scattering 
response  makes  the  standoff  detection  of  ground  targets  a  challenge,  especially  when  the  targets 
are  non-metallic  and  constraints  may  be  imposed  on  the  maximum  power  that  the  radar  transmitter 
can  radiate. 
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Figure  1.  Simulation  domain:  forward-looking  sensing  of  a  terrain  scene.  Note  that  the  vehicle,  transmitting  and 
receiving  antennas  are  not  included  in  the  finite -difference  time  domain  (FDTD)  grid — the  transmitting 
antennas  are  replaced  with  equivalent  aperture  current  distributions,  and  the  electric  fields  are  simply 
sampled  at  the  locations  of  the  receiving  antennas  for  image  generation. 

In  general,  as  shown  by  Liao  (201 1)  in  a  study  focused  mainly  on  landmines  and  metallic  shells, 
the  radar  responses  from  shallow-buried  targets  are — depending  on  the  polarization — at  least  8  to 
20  dB  weaker  than  those  from  their  on-surface  counterparts  and,  hence,  are  much  more  vulnerable 
to  interference  scattering  effects  arising  from  interface  roughness  and  subsurface  clutter.  As 
such,  the  understanding  of  the  scattering  responses  of  distributed  and  discrete  ground  clutter  (e.g., 
surface  unevenness,  soil  inhomogeneities,  plant  roots,  and  ground  debris)  is  expected  to  play  an 
important  role  in  assessing  the  detectability  limit  of  the  radar  for  sensing  buried  targets.  The 
scattering  interactions  within  a  radar  scene,  inherently,  are  intricate  functions  of  the  radar 
geometry  (i.e.,  transmitter/receiver-to-target  distance,  array  configuration),  transmission  and 
reception  polarizations,  frequency,  as  well  as  the  electrical  and  physical  properties  of  the  targets 
and  the  terrain  environment.  From  the  electromagnetic  perspective,  in  analyzing  a  low-to-ground 
system  such  as  the  SIRE  radar,  as  the  propagation  paths  defining  the  radiowave  interactions 
between  the  radar  transceiver  and  the  targets  adhere  to  the  ground,  existing  analytical  and  heuristic 
ground- surface  and  target-scattering  models  often  prove  to  be  insufficient  in  characterizing 
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various  backscattering  phenomena  of  interest.  Accordingly,  from  the  system  development 
perspective,  full-wave  electromagnetic  simulations  are  desired — and  necessary — for  an  accurate 
prediction  of  the  performance  of  the  forward-looking  radar. 

Previous  efforts  in  systematically  modeling  the  operation  of  radio  sensing  technology  for  the 
detection  of  ground  targets  include  the  studies  outlined  by  Bourgeois  and  Smith  (1996);  Bourgeois 
and  Smith  (1998);  Teixeira  et  al.  (1998);  Giirel  and  Oguz  (2000);  Giirel  and  Oguz  (2001);  Lampe 
and  Hooliger  (2001);  and  Uduwawala  et  al.  (2005) — albeit  the  focus  of  these  works  is  mostly  on 
the  downward-looking  sensing  geometry,  and  thus  the  simulation  domains  considered  therein  are 
relatively  small.  In  the  work  by  Teixeira  et  al.  (1998),  a  FDTD  method  is  used  for  computing  the 
scattering  from  metallic  and  plastic  pipes  buried  in  a  ground  with  dispersive  soil.  To  mitigate  the 
direct  signal  coupling  between  the  transmitter  and  receiver,  in  the  work  by  Giirel  and  Oguz  (2000), 
the  transmitting  and  receiving  elements  are  housed  in  a  shielding  structure  coated  with  absorbers; 
surface  roughness  and  soil  heterogeneities  are  also  included  in  the  model  by  populating  random 
scatterers  within  the  FDTD  computational  grid.  In  both  studies  above,  the  transmitters  and 
receivers  are  assumed  to  be  electric  dipoles;  UWB  radiators  more  akin  to  those  used  in  actual  GPR 
systems  are  analyzed  by  Bourgeois  and  Smith  (1996),  Lampe  and  Hooliger  (2001),  and 
Uduwawala  et  al.  (2005),  where  bow-tie  antennas  are  exploited  as  sources  of  illumination.  In  all 
the  aforementioned  studies,  the  emphasis  is  on  characterizing  the  ground  and  target-scattering 
responses,  and  therefore  no  imaging  results  are  included.  Note  that  in  addition  to  FDTD 
techniques,  many  method-of-moments  (MoM)-based  approaches  also  have  been  proposed  (He  et 
al.,  2000;  El-Shenawee  et  al.,  2001a;  El-Shenawee  et  al.,  2001b;  Johnson  and  Burkholder,  2001; 
El-Shenawee  and  Rappaport,  2002;  Johnson,  2002;  Wang  et  al.,  2003;  Altuncu  et  al.,  2006;  Guan 
et  al.,  2009)  to  facilitate  GPR-related  numerical  investigations — including  cases  where  the 
ground- surface  roughness  needs  to  be  taken  into  account;  sensing  and  imaging  at  grazing-angles, 
however,  are  not  discussed  in  detail  within  these  works. 

An  effort  by  the  authors  to  characterize  the  detection  operation  of  the  SIRE  radar  using  a  full-wave 
electromagnetic  solver  is  documented  in  previous  works  (Liao  et  al.,  2010;  Liao  and  Dogaru, 
2011a;  Liao  and  Dogaru,  2011b;  Liao  and  Dogaru,  2012),  in  which  only  far-field  interactions  are 
studied — that  is,  the  excitation  of  the  target  and  the  environment  is  supplied  by  pure  plane  waves, 
and  the  receivers  of  the  scattered  signal  are  placed  at  infinity.  In  the  current  work,  to  achieve  a 
more  accurate  emulation  of  the  sensing  operation  employed  by  the  SIRE  radar,  a  “near-field” 
LDTD  approach  is  proposed  to  provide  system-level  analysis,  and  for  tracking  the  propagation  and 
scattering  of  the  radio-wave  energy  through  the  computational  domain.  This  study  takes 
advantage  of  the  computing  resources  received  through  the  three-year  High  Performance 
Computing  Modernization  Program  Challenge  Project  award  entitled  “Advanced  Computational 
Modeling  of  Lorward-Looking  Radar  for  the  Detection  of  Landmines  and  IEDs.”  This  report 
provides  a  summary  of  the  works  performed  during  the  first  year  of  the  project  and  is  organized  as 
follows:  section  2,  the  framework  of  the  full-wave-based  radar  emulator  is  detailed;  section  3,  an 
investigation  of  the  validity  of  the  “near-field”  FDTD  algorithm  is  put  forth;  and  section  4,  a 
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demonstration  of  the  effects  of  ground  clutter  on  target  imaging  generated  by  the  time-reversal 
technique  is  presented  with  the  simulation  of  a  large  terrain  scene  consisting  of  targets  in  a  rough 
ground  environment.  In  sum,  this  study  illustrates  the  practicality  of  the  chosen  full-wave 
simulation  strategy  for  the  emulation  of  forward-looking  radar  operation  and  imaging. 


2.  Electromagnetic  Simulation  and  Analysis  Framework 


In  demonstrating  the  modeling  framework  put  forth  in  this  study,  electromagnetic  simulation  and 
coherent  imaging  of  the  scene  depicted  in  figure  1  is  considered.  Note  that  during  in-field 
operations,  the  ARL  SIRE  radar  integrates  scattering  responses  gathered  at  multiple  vehicle 
locations  to  form  the  final  image  for  a  scene.  As  this  is  the  first  phase  of  the  numerical 
investigation,  in  order  to  establish  a  baseline  for  validation  and  identify  any  potential  fundamental 
limitations  of  the  solver,  the  radar  transceiver  is  assumed  to  be  fixed  in  position,  with  the  image 
synthesized  from  return  signals  captured  at  only  a  single  aperture.  With  these  considerations,  the 
components  of  the  forward-looking  radar  emulation  engine  are  outlined  below. 

2.1  Terrain  Characterization 

The  ground  is  modeled  as  a  non-dispersive,  non-permeable  homogeneous  medium  with  effective 
relative  dielectric  constant  er  and  conductivity  a,/.  The  electrical  properties  of  the  ground  are 
dependent  upon  the  soil  composition  and  moisture  content;  er  and  a,i  can  be  either  determined 
using  direct  dielectric  measurement  techniques  or  approximated  with  empirical-based 
formulations  (Dobson  et  al.,  1985;  Peplilnski  et  al.,  1995).  The  surface  statistics  of  the  ground,  in 
general,  are  described  by  two  functions:  the  probability  density  function  of  the  height  variations 
and  the  surface  autocorrelation  function.  Within  the  scope  of  this  study,  a  two-dimensional 
zero-mean  surface  profile  obeying  Gaussian  statistics  {hrms,  lc} — where  hrms  is  the  root  mean 
square  (rms)  height  and  lc  the  correlation  length — is  assumed  and  generated  from  its  randomized 
power  spectral  density  function  by  following  the  procedure  prescribed  by  Thorsos  (1988). 

2.2  Ground  Targets 

As  illustrated  in  figure  1,  targets  in  the  form  of  landmines  (either  anti-tank  or  anti-personnel, 
metallic  or  plastic,  buried  or  on-surface)  and  155-mm  shells  (metallic,  buried)  are  considered. 
Specifically,  these  targets — the  shapes,  dimensions,  orientations,  and  locations  of  which  are  also 
indicated  by  the  white  outlines  in  figure  4a — include  the  following:  (1)  Buried  metallic 
anti-personnel  landmine;  (2)  buried  plastic  (£y=3.1,  oj=2  mS/m, )  anti-personnel  landmine; 

(3)  on-surface  plastic  (c,=3. 1 ,  at{=2  mS/m )  anti-personnel  landmine;  (4)  buried  metallic  155-mm 
shell;  (5)  buried  metallic  anti-tank  landmine;  (6)  on-surface  metallic  anti-tank  landmine;  (7)  buried 
metallic  155-mm  shell;  (8)  Buried  metallic  155-mm  shell;  (9)  On-surface  metallic  anti-personnel 
landmine;  (10)  buried  plastic  (£,=3.1,  Od=  2  mS/m )  anti-tank  landmine;  and  (11)  on-surface  plastic 
(8,= 3.1,  (7,1=2  mS/m )  anti-tank  landmine.  Buried  targets  are  positioned  at  3  cm  below  the  surface. 
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2.3  Transmitting  and  Receiving  Antenna  Treatment 

The  SIRE  radar  employs  two  transverse  electromagnetic  (TEM)  horn  antennas  as  transmitters  and 
16  Vivaldi  antennas  distributed  over  a  2  m-wide  aperture  as  receivers.  In  the  current 
implementation  of  the  algorithm,  the  TEM  hom  transmitting  antennas  are  simply  replaced  by 
equivalent  current  distributions;  that  is,  the  following  electric  and  magnetic  current  sheets  are 
positioned  at  the  horn  apertures  (figure  1):  J  =  xxH\  M  =  Exx  ,  where  the  electric  and  magnetic 
fields  are  assumed  to  be  constant  over  the  apertures.  The  Vivaldi  antennas  are  not  included  in  the 
FDTD  model:  the  electric  fields  are  sampled  at  the  locations  of  the  16  antennas  for  image 
generation. 

2.4  Three-dimensional  “Near-Field”  FDTD  Algorithm 

The  composite  scattering  from  the  targets  and  the  rough  ground  is  calculated  numerically  with  the 
FDTD  method  (Taflove  and  Hagness,  2005),  which  directly  solves  the  differential  forms  of 
Maxwell’s  equations  in  the  spatial  and  temporal  domains.  In  essence,  the  three-dimensional 
computational  domain  is  decomposed  into  small  cubic  cells,  with  each  cell  assigned  the  material 
dielectric  properties  of  the  corresponding  region  in  space.  Electric  and  magnetic  field 
components  are  sampled  along  the  cells’  edges  and  faces,  and  the  field  spatial  distribution  is  then 
propagated  in  time  according  to  a  time-marching  scheme  that  implements  a  discretized  version  of 
Maxwell’s  equations.  The  entire  computational  grid  is  terminated  by  a  perfectly  matched  layer 
(PML)  that  emulates  reflection-less  propagation  of  the  outgoing  waves  in  infinite  space.  In  the 
current  application,  an  infinite  half-space  (air/ground  layered  media)  is  considered  as  the 
background  environment  in  which  the  targets  are  embedded;  the  rough  ground- surface  is 
introduced  as  a  finite,  random  perturbation  of  the  flat  interface  between  the  two  media  (the 
approach  is  similar  to  that  described  by  Liao  and  Dogaru  [2012]).  The  overall  FDTD  code  has 
been  fully  parallelized  and  runs  on  large  distributed  computer  systems  using  the  Message  Passing 
Interface  (MPI)  framework:  the  computational  domain  is  decomposed  into  rectangular 
sub-domains  and  the  FDTD  equations  are  solved  separately  for  each  sub-domain  within  one  MPI 
process.  Since  only  the  electric  and  magnetic  field  samples  in  the  boundary-layer  between  two 
adjacent  sub-domains  need  to  be  exchanged  between  the  respective  processes,  the  code  is  highly 
scalable,  even  when  the  number  of  MPI  processes  is  in  the  hundreds. 

The  field  excitation  sources  in  the  algorithm  are  implemented  as  electric  and  magnetic  current 
distributions  in  space  (“soft  sources”),  which  are  added  to  the  regular  FDTD  field  update  equations 
at  the  desired  sample  points.  Note  that  by  choosing  this  approach,  direct,  full-scale  modeling  of 
the  transmitting  antenna  inside  the  FDTD  code  is  avoided;  instead,  the  antenna  is  replaced  with 
electric  and/or  magnetic  current  distributions  by  employing  the  equivalence  principle. 
Computation  of  these  current  distributions  is  performed  externally  and  can  be  implemented  in 
various  ways — based  on  approximate  analytical  formulations  or  exact  numerical  analysis.  The 
frequency  domain  equivalent  currents  replacing  the  transmitting  antenna  are  subsequently 
modulated  by  the  spectrum  of  an  UWB  pulse  (in  this  case,  a  fourth-order  Rayleigh  pulse),  and  their 
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time-domain  versions  are  used  as  excitation  sources  in  the  FDTD  time-marching  scheme.  In  the 
numerical  examples  presented  in  sections  3  and  4,  the  sources  are  implemented  either  as  an 
infinitesimal  electric  dipole,  or  as  uniform  electric  and  magnetic  current  sheets  distributed  on  a 
rectangular  aperture  (the  latter  approximate  the  currents  produced  by  a  TEM  horn  antenna,  as 
explained  in  section  2.3).  During  the  simulation  run,  the  electric  fields  at  the  receiving  antenna 
locations  are  recorded  for  image  generation;  specifically,  the  frequency  domain  version  of  the 
received  electric  field  components  is  obtained  by  deconvolving  the  excitation  pulse  from  the 
FDTD  simulation  results. 

For  sensing  applications  at  low-depression-angles,  one  salient  feature  of  field  interactions  in  the 
vicinity  of  the  air/ground  interface  is  the  near  cancellation  of  the  direct  and  reflected  waves — as  it 
occurs  during  the  illumination  of  the  target  by  the  transmitter,  as  well  as  during  the  propagation  of 
the  backscattered  target  signal  to  the  receiver.  Consequently,  the  direct  fields  that  propagate 
between  closely  spaced  transmitting  and  receiving  points  can  be  many  orders -of-magnitude  larger 
than  those  received  via  scattering  from  the  target.  As  such,  to  properly  characterize  the  complete 
field  behavior,  stringent  requirements  must  be  placed  on  the  precision  and  the  dynamic  range  of 
the  solver  to  ensure  very  weak  scattering  responses  are  captured.  Given  these  considerations,  the 
only  solution  for  obtaining  accurate  numerical  values  of  the  fields  scattered  by  a  target  with  the 
“near-field”  code  is  to  run  the  FDTD  time-marching  process  twice — once  for  the  configuration 
that  includes  the  target  and  once  without  the  target — and  then  subtract  the  received  fields  from  the 
two  cases.  The  final  result  is  the  “perturbation”  in  the  total  field  introduced  by  the  presence  of  the 
target,  which  is  exactly  the  definition  of  the  scattered  field.  Note  that  some  types  of  errors 
associated  with  the  direct  (or  “unperturbed”)  field  are  completely  cancelled  by  this  procedure. 
However,  other  errors,  such  as  numerical  dispersion  artifacts,  cannot  be  eliminated  this  way  and 
must  be  kept  under  control  by  choosing  a  sufficiently-fine  grid  discretization. 

2.5  Time-reversal  Imaging  With  Asymptotics 

A  coherent  field  integration  technique  is  employed  to  process  the  FDTD-computed  scattered 
signals  to  obtain  an  image  of  the  illuminated  scene.  According  to  the  time-reversal,  or 
phase-conjugation,  method,  since  the  received  scattered  field  can  be  written  as 

Es={fR,fT,o})  =  Y{o})G(f5,fT,co)G{fR,rs,a>),  ^ 

an  approximate  objective  (image)  function  can  be  formed  with  the  following  expression  (Borcea  et 
al.,  2003;  Sarabandi  et  al.,  2004): 

0  =  {r)  =  YL^E*s(fR,rT,(o)G(f,rT,co)G{fR,r,(o)\ 

where  the  asterisk  notation  denotes  the  phase  conjugation  operation;  G(f,T',co)  is  the  Green’s 

function  of  the  environment  ( r  is  the  position  of  the  observation  point;  f  is  the  position  of  the 
source  point);  rT,rR,  and  r  are  the  locations  of  the  transmitter,  receiver,  and  scatterer,  respectively; 

Y(co)  is  the  “radiation  admittance”  of  the  scatterer — a  measure  of  how  easily  currents  can  be 
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induced  throughout  the  object.  The  cross-  and  down-range  resolution  of  the  image  are 
determined  by  the  sensor  aperture  and  system  bandwidth,  respectively,  of  the  coherent  summation 
performed  in  equation  2.  For  sensing  in  the  presence  of  a  randomly- varying  ground  interface,  as 
the  exact  propagation  Green’s  function  is  not  known,  G(f,f’,a> )  is  approximated  by  the 

half-space  Green’s  function;  as  such,  the  target  image  is  corrupted  by  ground  clutter. 

For  the  computation  of  G  (r ,  co)  — or  for  the  characterization  of  radiation  from  current  sources  in 
the  presence  of  a  half-space,  note  that  exact  closed-form  expressions  can  be  formulated  in  terms  of 
integrals  of  the  Sommerfeld  type;  however,  these  expressions  are  difficult  and  time  consuming  to 
evaluate  numerically  (especially  when  the  source  and  observation  points  are  far  apart).  To 
completely  avoid  the  use  of  Sommerfeld  integrals,  the  transformation  of  these  integrals  into 
asymptotic  forms  can  be  carried  out  by  applying  the  method  of  steepest  descents  (Liao  and 
Sarabandi,  2005).  The  generalized  asymptotic  expressions — simplified  to  the  half-space 
propagation  scenario — from  the  work  of  Liao  (2009)  are  employed  here  in  equation  2  to  achieve  a 
more  efficient  imaging  process. 


3.  Validation  of  FDTD  Solver 


The  accuracy  of  the  simulation  approach  outlined  in  section  2  is  determined  by  comparing 
FDTD-derived  results  with  those  obtained  using  another  full-wave  technique,  namely,  MoM,  as 
implemented  in  the  commercially  available  electromagnetic  software  FEKO.  Only  relatively 
simple  canonical  problems  are  considered  in  the  validation  process  here,  since  more  complicated 
scenarios  become  intractable  in  FEKO. 

The  images  generated  from  FDTD-  and  MoM-computed  scattering  responses  for  an  on-surface 
target  are  compared  in  figures  2  and  3.  A  small  perfectly  conducting  cube  (side  length=88  mm)  is 
considered  in  figure  2,  while  a  perfectly  conducting  cylinder  (radius=height=15  cm)  is  considered 
in  figure  3.  The  antenna  source  is  an  infinitesimal  electric  dipole  (vertically-polarized)  located 
above  the  ground  at  2  m  height;  16  evenly  spaced  receiving  points  are  distributed  over  a  2-m-wide 
aperture  parallel  to  the  y-axis  with  the  receiving  aperture  centered  immediately  below  the  dipole,  at 
a  height  of  1.9  m.  The  targets  are  positioned  at  a  down-range  distance  of  either  3.24  m  or  18.24  m 
from  the  source. 
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Figure  2.  Comparison  of  FDTD-  and  MoM-simulated  imaging  results  for  an  on-surface,  perfectly-conducting  cube 
at  two  ranges:  (a)  FDTD:  range=3.24  m;  (b)  MoM:  range=3.24  m;  (c)  FDTD:  range=18.24  m;  and 
(d)  MoM:  range=  18.24  m.  Ground  properties:  e,=6,  <t„=  1 0  mSlm.  Frequency  span:  0.1-1. 4  GHz. 
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Figure  3.  Comparison  of  FDTD-  and  MoM-simulated  imaging  results  for  an  on-surface,  perfectly-conducting 

cylinder  at  two  ranges:  (a)  FDTD:  range=3.24  m;  (b)  MoM:  range=3.24  m;  (c)  FDTD:  range=18.24  m; 
and  (d)  MoM:  range=18.24  m.  Ground  properties:  a,j=  1 0  mS/m.  Frequency  span:  0. 1-1.4  GHz. 

As  conveyed  in  the  figures,  the  FDTD  results  are  in  good  agreement  with  those  from  MoM  in 
characterizing  the  primary  return  from  the  target  for  all  test  cases  analyzed.  Unfortunately,  the 
FDTD  response  contains  a  late-time  return  (secondary  “echo”)  that  does  not  appear  in  the  MoM 
solution.  It  is  conjectured  here  that  this  late -time  “echo,”  which  lies  23  to  30  dB  below  the  main 
response,  is  purely  an  FDTD  simulation  artifact — a  byproduct  of  numerical  dispersion  inherent  in 
FDTD  and  of  the  imperfectness  of  the  absorbing  boundary  condition.  It  is  also  observed  that  the 
magnitude  of  this  error  can  be  reduced  by  decreasing  the  grid  cell  size.  This  has  to  be  done,  of 
course,  at  the  expense  of  increasing  the  central  processing  unit  (CPU)  and  memory  usage. 
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4.  Large-scale  Simulation  and  Imaging  Results 


Full -wave  characterization  of  the  scene  shown  in  figure  1  is  performed  with  the  “near-field” 
FDTD  algorithm  described  in  section  2.  The  computational  domain  has  dimensions  of 
20  m  x  10mx4m  and  is  discretized  into  2.6  billion  FDTD  cells.  Parallelized  simulations  are 
performed  after  partitioning  the  scene  into  128  sub-domains  (i.e.,  a  total  of  128  processors  are 
employed  in  the  MPI  network).  Approximately  1,000  CPU-hours  are  needed  for  each  simulation 
run.  As  mentioned  in  section  2.4,  the  computation  of  the  scattered  field  essentially  requires  two 
runs:  one  for  the  “unperturbed”  problem — the  source  in  the  presence  of  the  flat  ground  only, 
without  targets — and  one  for  the  “perturbed”  problem — with  targets  and  any  ground- surface 
roughness.  In  view  of  the  fact  that  the  SIRE  system  employs  two  transmitters,  two  separate  sets 
of  scattering  responses  (corresponding  to  two  sources  of  excitation)  are  synthesized  together  to 
form  the  image  for  the  scene — that  is,  the  transmitters  are  excited  sequentially,  not  simultaneously, 
for  each  simulation. 

The  effects  of  a  rough  ground- surface  on  focusing  are  demonstrated  with  the  simulation  results 
displayed  in  figures  4  and  5,  which  show  the  co-polarization  images  of  the  scene  for  different 
degrees  of  surface  roughness:  (a)  flat  ground;  (b)  hrms=  1.2  cm,  4=14.93  cm;  and  (c)  hrms=  1-6  cm, 
4=14.93  cm.  It  should  be  emphasized  that  the  incident  (direct)  field  amplitude  at  the  location  of 
each  image  pixel,  as  well  as  the  propagation  distance  from  the  pixel  to  the  receiver  have  been  taken 
into  account  in  the  imaging  process.  In  effect,  each  derived  image  can  be  interpreted  to  be  a  map 
containing  the  spatial  variation  of  the  average  bi-static  radar  cross  section  (RCS)  a(rs)  for  the 

scene,  in  which  the  averaging  operation  is  carried  over  multiple  transmitter/receiver  locations  and 
a  range  of  frequencies.  Consequently,  consistent  with  this  RCS  interpretation  of  the  observed 
image  intensities,  a  particular  scatterer  would  appear  brighter  the  closer  it  is  to  the  radar — due  to 
illumination  at  a  larger  depression  angle  (i.e.,  away  from  grazing-angle)  that  occurs  in  the 
near-range.  (Depending  on  the  “normalization”  method  used,  the  image  can  also  be  construed  as 
a  map  of  the  electric  surface  and  volumetric  current  densities  induced  throughout  the  scene  by  the 
source;  this  is  realized  by  noting  that  the  term  Y ( co)G(rs,rT,co )  in  equation  1  is  basically  a 

transformation  of  the  incident  fields  into  currents.)  From  the  analysis  of  the  point-spread  function 
(see  the  appendix),  note  that  the  sensing  resolution  for  the  current  setup  is  c/2A/in  the  down-range 
direction  and  XcrlLeff  in  the  cross-range  direction,  where  c  denotes  the  speed  of  light;  A/ is  the 
frequency  bandwidth;  Xc  is  the  wavelength  at  the  center  frequency;  r  is  the  range  distance;  and  Leff 
is  the  effective  aperture  length — in  this  case,  it  is  the  sum  of  the  lengths  of  the  transmitting  and 
receiving  apertures. 
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Figure  4.  Simulated  vertical-vertical  (vv)-polarized  images  for 
20  mxlO  m  scene:  (a)  Ground  with  flat  surface; 

(b)  Ground  with  randomly  rough  surface,  1 .2  cm, 
lc=  14.93  cm;  (c)  Ground  with  randomly  rough  surface, 
hnm=  1.6  cm,  l,=  1 4.93  cm.  Ground  electrical  properties: 
e,=6,  r>,/=  1 0  mS/m.  Frequency  span:  0.3-1. 5  GHz. 
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Figure  5.  Simulated  horizontal-horizontal  (lili)- polarized  images 
for  20  mxlO  m  scene:  (a)  Ground  with  flat  surface; 

(b)  Ground  with  randomly  rough  surface,  hrms  =1.2  cm, 

/,= 14.93  cm;  (c)  Ground  with  randomly  rough  surface, 
hrms=  1 .6  cm,  lc=  14.93  cm.  Ground  electrical  properties: 
e,=6,  <T,i=  1 0  mS/m.  Frequency  span:  0.3-1. 5  GHz. 

As  seen  in  the  figures,  buried  targets  are  generally  more  difficult  to  discern  from  the  background  as 
compared  to  on-surface  targets.  The  response  from  the  larger  targets  (anti-tank  landmines)  in  the 
backscattering  directions  is  primarily  due  to  edge-scattering  from  the  front-  and  back-rims  of  the 
top  of  the  cylindrical  structure;  the  scattering  centers  corresponding  to  these  two  edges  are 
resolved  in  the  images — for  both  on-surface  and  buried  cases.  For  the  smaller  targets 
(anti-personnel  landmines),  the  two  scattering-edges  tend  to  merge  into  one,  since  the  imaging 
configuration  as  chosen  here  is  unable  to  provide  sufficient  resolution  to  separate  the  two  phase 
centers.  In  the  presence  of  a  rough  ground  interface,  these  scattering  effects  are  less  apparent  in 
the  images  for  the  buried  targets  due  to  interference  from  surface  scattering.  A  buried  plastic 
target  (either  large  or  small  landmine)  is  especially  hard  to  detect,  owing  to  the  limited  dielectric 
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contrast  between  the  target  and  the  soil  background.  (The  FDTD  artifact  discussed  in  section  3 
can  also  be  seen  in  these  images,  when  its  intensity  falls  within  the  dynamic  range  shown — which 
is  true  for  targets  5  and  6  at  vv,  for  example.)  The  backscattering  return  from  the  ground-surface 
itself  is  stronger  at  vv  than  hh,  and  the  response  at  both  polarizations  tends  to  decrease  as  the 
incidence  angle  approaches  grazing.  These  observations  are  consistent  with  previously  reported 
results  derived  from  far-field  simulations  (Liao  and  Dogaru,  2011b)  and  scatterometer 
measurements  (Oh  et  al.,  1992). 

Comparing  all  target  responses  across  the  two  polarizations,  on  average,  vv  provides  the  stronger 
return,  and  it  exceeds  hh  by  about  5  to  10  dB.  As  a  vertically  polarized  signal  experiences  less 
transmission  loss  at  the  air-ground  interface  and  suffers  less  propagation  loss  in  general,  the  above 
observation  supports  the  notion  that  vv  is  the  more  desirable  configuration  for  detecting  targets 
either  on  top  of  or  buried  under  a  flat  ground — at  the  low-frequency  region  studied  (<1.5  GHz). 

In  the  presence  of  ground-surface  roughness,  detection  performance  must  be  evaluated  with  a 
consideration  of  the  target-to-background  clutter  ratio.  As  such,  based  on  the  results  of  these 
figures,  given  that  the  surface  response  appears  as  a  stronger  interference  component  at  vv,  it  is  not 
apparent  that  sensing  with  vertically-polarized  waves  provides  a  definitive  advantage — for  the  set 
of  parameters  examined  herein.  The  statistical  distribution  of  the  ground  polarimetric  scattering 
response  as  a  function  of  the  surface  properties  can  also  be  estimated  to  facilitate  a  more 
comprehensive  analysis  of  the  target  detection  performance. 


5.  Conclusion 


To  predict  the  performance  of  forward-looking  radar  technology  deployed  in  complex 
environments,  large-scale,  full-wave  analysis  of  multi-static  target-imaging  in  a  rough-ground 
environment  is  described.  The  emulation  methodology  employs  a  parallelized  three-dimensional 
“near-field”  FDTD  algorithm  in  characterizing  the  electromagnetic  scattering  from  the 
ground- surface  and  buried  and  on-surface  targets.  The  transmitting  antennas  are  represented  by 
equivalent  current  sheets  placed  at  the  real  antenna  aperture  locations.  Subsequent  focusing  of 
the  scattered  fields  at  the  receiver  locations  into  an  image  is  obtained  by  exploiting  the 
time -reversal  technique,  with  the  aid  of  asymptotics.  The  overall  simulator  is  expected  to  be  a 
useful  tool  for  studying  the  effects  of  ground  properties — both  physical  and  electrical — on  the 
signatures  of  various  targets. 

As  evident  from  the  simulation-based  validation  studies,  the  primary  scattering  returns  of  the 
targets  are  accurately  characterized  by  the  proposed  FDTD  approach;  however,  numerical 
artifacts — which  manifest  as  “echoes”  of  the  main  imaging  response — are  also  observed.  Future 
work  includes  investigating  the  fundamental  causes  of  this  error  component  as  well  as  exploring 
means  to  mitigate  its  effects  without  compromising  the  computational  efficiency  of  the  solver. 
Efforts  in  extending  the  algorithm  to  accommodate  the  full  models  of  the  transmitting  and 
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receiving  antennas  will  also  be  undertaken.  In  addition,  multi-aperture  scattering  and  imaging 
simulations  will  be  performed  in  the  next  phase  of  the  study  to  provide  a  true  emulation  of  the 
forward-motion-based  operational  capability  of  the  SIRE  system. 
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Appendix.  Derivation  of  the  Imaging  Point-Spread  Function 


Consider  an  array  of  vertically  polarized  point  radiators  and  a  point- sc atterer  directly  on  the 
ground  at  rs  =(rs,9s,<j>s),  approximating  the  Green’s  function  for  the  interactions  of  the  co-polarized 

field  component  as  (Liao,  2009) 


and  assuming  observation  with  an  array  of  isotropic  point  receivers,  the  point-spread  function 
(PSF)  can  be  constructed  from  equation  2  as 


PSF(r,  rs )  «  'F  (G  0"*’ f. lp  {?,  fT)e~ 


(A-2) 


Using  far- field  expressions  for  the  distance  vectors,  the  F  (■)  terms  can  be  moved  outside  the 
summation  operation,  and  orienting  the  transmitting  and  receiving  arrays  parallel  to  the  y-axis 
(xr=.xt=0),  the  function  in  equation  A-2 — on  the  xy-plane — becomes 


PSF  ( r ,  rs )  cc  F  ( rs ,  rT )  F  ( rR ,  rs )  F  ( r ,  rT )  F  ( rR ,  r ) . 

y  y  y  ^2 ink8k{rs  —  r )nRSyR (sin ^-sin </>s )^/(nit^fe-i-fec)n7-^yT(sin^— sin^s ) 


(A- 3) 


where  5k ,  dyR,  and  civ/  are  the  k- space  sampling  interval,  the  receiver  spacing,  and  the  transmitter 
spacing,  respectively;  kc  is  the  center  wavenumber;  the  wavenumber  bandwidth  is  A k=Nu5k\  and 
the  aperture  sizes  are  Lr=NrSvr,  and  L-i-NtSvt-  Applying  a  stationary-phase-like  approach  to 
evaluate  the  summation  in  k-space  and  noting  the  relation 


y  g-im  =  sin(f) 

n'  sin(t)  ’ 


(A-4) 


the  amplitude  of  equation  A-3  approximately  simplifies  to 

|P5F(r,rs)|oc|F(rs,4)F(rfl,rs)F(r,rr)F(rfl,F)|. 

sin(M(r,-r))  ^  sln(  _  /  L,-k,( sin^-sin^,) 

sin(«(r,-r))  j  (  2 


(A- 5) 


where  it  is  assumed  that  there  is  a  transmitter  at  each  end  of  the  aperture  Lj.  The  first  sin  (• )  /  sin 
(• )  term  in  equation  A-5  dictates  the  behavior  of  the  PSF  in  range:  it  can  be  shown  that  the  first  null 
occurs  at  c/2Af- — which  is  the  resolution  in  range,  dr,  independent  of  rs  and  <f>s.  The  periodicity  of 
this  term  determines  the  frequency  sampling  interval  needed:  for  maximum-range  rsjnax,  the 
condition  5f<c/2rs  max  must  be  met.  Similarly,  the  cross-range  properties  of  the  PSF  are  controlled 
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by  the  last  two  terms  in  equation  A-5.  For  the  current  two-transmitter  system  (and  L=LT=LR), 
based  on  the  location  of  the  first  null  of  these  two  terms,  it  can  be  readily  seen  that  the  angular 
resolution  depends  on  the  angular  position  of  the  scatterer  as 


8(j)  =  (j)s  -  sin 


°^<§- 


(A-6) 


At  (f>s= 0,  S(f>~/.,J2L:  the  angular  resolution  degrades  as  the  scatterer  moves  away  from  the  broadside 
direction  of  the  aperture.  In  general,  the  spatial  sampling  should  be  less  than  2min/2;  however,  for 
an  UWB  system,  it  has  been  shown  (Zhuge  and  Yarovoy,  2011)  that  high  cross-range  resolution 
can  be  achieved  with  coarser  sampling  than  the  stated  criterion. 
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List  of  Symbols,  Abbreviations,  and  Acronyms 


ARL 

U.S.  Army  Research  Laboratory 

CPU 

central  processing  unit 

FDTD 

finite-difference  time-domain 

GPR 

ground-penetrating  radar 

hh 

horizontal-horizontal 

IED 

improvised  explosive  devices 

MoM 

method-of-moments 

MPI 

message-passing  interface 

PML 

perfectly-matched  layer 

PSF 

point-spread  function 

RCS 

radar  cross  section 

rms 

root  mean  square 

SIRE 

synchronous  impulse  reconstruction 

TEM 

transverse  electromagnetic 

UWB 

ultra-wideband 

vv 

vertical-vertical 
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